cd "L:/Project-Fire_Perimeter/Fire_perimeter/Processing/Ignitions/BufferSections/Sections/"

capture program drop mymargins
program mymargins

{
	syntax varlist, estname(string)
	local ll: display %9.0fc `e(ll)'	
	margins, dydx(`varlist') post
	est sto `estname'
	estadd local ll "`ll'"
}	
end

capture program drop distlag
program distlag

syntax namelist [if] , method(string) params(numlist) tmethod(string)  ///
	[recode BURNTmeasure(string) tcuts(numlist) lagvars(namelist) costlags(namelist) ///
	save(string) fe(namelist) cluster(namelist) output(string)] 
			
	qui {
	
	/************************************************************************************/
	/*Recode burnt variable. Code as zero all sections that are burnt after a section  	*/
	/*is unburnt. This guarantees for each fire direction will have a single burn 		*/
	/*duration, is an imperfect way of dealing with convex fires.						*/
		
	if "`burntmeasure'" == "centroid" {
		drop burnt
		rename burntalt burnt
		}
	else drop burntalt
		
	if "`recode'" != "recode" { /*Code as zero if beyond first extinguished cell*/
		sort firenum directionnum distnum
		capture drop count
		bysort firenum directionnum: gen count = sum(burnt)
		replace burnt = 0 if count < distnum & burnt == 1
		drop count
		}
	else if "`recode'" == "recode" { /*Code as 1 if not last extinguished cell*/
		sort firenum directionnum distnum
		capture drop count
		bysort firenum directionnum: gen count = sum(burnt)
		bysort firenum directionnum: egen max = max(count)
		replace burnt = 1 if count < max
		drop count max
		}
	
	sort firenum directionnum distnum	
	bysort firenum directionnum: egen x = sum(burnt)
	bysort  firenum directionnum: egen duration = max(x)
	drop x
		
	gen y = (distnum == duration + 1)
	
	/**********************************************************************************/
	/* Generate binned time-difference variables									  */		
	

	rename t toa
	
	if "`tmethod'" == "bins" {

		local q = 0
		if "`tcuts'" != "" {
			foreach cut in `tcuts' {
				local q = `q' + 1
				if `q' == 1 gen t = `cut' if lnt <= `cut'
				else replace t = `cut' if lnt > `tmin' & lnt <= `cut'
				local tmin = `cut'
				}
			tab t, gen(t)
					
			local j = 0
			foreach var of varlist t1-t`q' {
				replace `var' = 0 if tmissing == 1
				local j = `j' + 1
				if `j' == 1 {
					local max: word `j' of `tcuts' 
					label var `var' "<`max'"
					}
				else if `j' > 1 {
					local k = `j' - 1
					local max: word `j' of `tcuts' 
					local min: word `k' of `tcuts'
					label var `var' "`min'-`max'"
					}
				}
			label var tmissing "Missing"
			}
		}
			
			
	/**********************************************************************************/
	/* Generate var list for unrestricted distributed lag models					  */
	
	
	if "`method'" == "unrestricted" {
		
		local q: word 1 of `params'
		
		sort firedirid distnum
		tsset firedirid distnum
		
		foreach var of varlist `lagvars' `costlags' {
			
			forvalues i = 0/`q' {
			
				gen F`i'`var' = F`i'.`var'
				local lags `lags' F`i'`var'
				
				}
			}	
		}
	 	
	/**********************************************************************************/
	/* Generate Almon variables														  */			

	if "`method'" == "almon" {
	
		local p: word 1 of `params'
		local n: word 2 of `params'
		
		sort firedirid distnum
		tsset firedirid distnum
	
		foreach var of varlist `lagvars' `costlags' {
						
			local cpos: list posof "`var'" in costlags
						
			forvalues i = 0/`p' {
				gen w`var'`i' = 0
				forvalues j = 1/`n' {
					replace w`var'`i' = cond(firedirid[_n] == firedirid[_n+`j'], w`var'`i' + `j'^`i'*`var'[_n + `j'], .)
					}			
				local lags `lags' w`var'`i'	
				}
			if "`cpos'" == "0" replace w`var'0 = w`var'0 + `var'
			}	
		}
	
	
	/**********************************************************************************/
	/* Generate linearly declining lag weights										  */			
	
	if "`method'" == "linear" {
	
		local q: word 1 of `params'
		
		sort firedirid distnum
		tsset firedirid distnum
		
		foreach var of varlist `lagvars' `costlags' {
		
			local cpos: list posof "`var'" in costlags
		
			gen w`var' = `var'
			
			forvalues s = 1/`q' {
				replace w`var' = cond(firedirid[_n] == firedirid[_n +`s'], w`var' + ((`q' + 1 - `s')/(`q' + 1))*`var'[_n + `s'], .)
				}
				
			local lags `lags' w`var'
			if "`cpos'" == "0" replace w`var' = w`var' + `var'
							
			}
		}
		
	/**************************************************************************************/
	/* Drop observations after fire has gone out										*/
	
	if "`dropout'" == "" | "`dropout'" == "dropout" {
		drop if distnum > duration + 1
		}
	else if "`dropout'" == "nodropout" {
		gen fireout = (distnum > duration + 1)
		}
		
	/**************************************************************************************/
	/* Run regressions																	*/	
		
	*Dependent variables
	local depvar: word 1 of `namelist'
	
	*Independent variables
	local indepvars `lags'
	if "`method'" != "unrestricted" local indepvars `lags' `costlags'
	local listsize: list sizeof namelist
	
	if `listsize' > 1 {
		forvalues i = 2/`listsize' {
			local var: word `i' of `namelist'
			local indepvars `indepvars' `var'
			}
		}
		
	di "`indepvars'"

	if "`fe'" == "" {
		noi reg `depvar' `indepvars' i.distnum `if', cluster(`cluster') 
		}
	if "`fe'" != "" {
		noi reghdfe `depvar' `indepvars' i.distnum `if', vce(cluster `cluster') absorb(`fe')
		}
				
	if "`output'" != "table" {
	
		if "`method'" == "unrestricted" {
			
			local word1 : word 1 of `lagvars'
			local q: word 1 of `params'
			
			foreach var of varlist `lagvars' `costlags' {

				*local cpos: list posof "`var'" in costlags

				forvalues j = 0/`q' { 
			
					lincom F`j'`var'
					
					if "`j'" == "0" & "`var'" == "`word1'" mat def blag = [r(estimate),r(se),r(lb),r(ub)]
					else mat def blag = [blag \ r(estimate),r(se),r(lb),r(ub)]

					}
				}
					
			local n = `q'
			
			}			
			
		if "`method'" == "linear" {
		
			local word1 : word 1 of `lagvars'
			local q: word 1 of `params'
		
			foreach var of varlist `lagvars' `costlags' {
			
				local cpos: list posof "`var'" in costlags

				forvalues j = 0/`q' { 
				
					if "`cpos'" != "0" & "`j'" == "0" lincom `var'

					else {
						local m = ((`q' + 1 - `j')/(`q' + 1))
						lincom `m'*_b[w`var']
						}
					
					if "`j'" == "0" & "`var'" == "`word1'" mat def blag = [r(estimate),r(se),r(lb),r(ub)]
					else mat def blag = [blag \ r(estimate),r(se),r(lb),r(ub)]
				
					}
				}
				
			local n = `q'
			
			}

		
		if "`method'" == "almon" {
			
			/*Calculate restricted lags and store in matrix*/
			
			local word1 : word 1 of `lagvars'
			local p: word 1 of `params'
			local n: word 2 of `params'
		
			foreach var of varlist `lagvars' `costlags' {
			
				local cpos: list posof "`var'" in costlags
					
				forvalues j = 0/`n' { 
					
					local b`j'`var'
					
					forvalues i = 0/`p' {
					
						local m = `j'^`i'
						if "`i'" != "`p'" local b`j'`var' `"`b`j'`var'' _b[w`var'`i']*`m' + "'
						else local b`j'`var' `"`b`j'`var'' _b[w`var'`i']*`m'"'
											
						}
						
					if "`cpos'" != "0" local b0`var' `var'
						
					di "`b`j'`var''"
					lincom `b`j'`var''
				
					if "`j'" == "0" & "`var'" == "`word1'" mat def blag = [r(estimate),r(se),r(lb),r(ub)]
					else mat def blag = [blag \ r(estimate),r(se),r(lb),r(ub)]
					
					}
				}
			}
				
			/* Create variables for figure */ 

			clear
			svmat blag
			rename blag1 beta
			rename blag2 se
			rename blag3 ul
			rename blag4 ll	
			gen var = ""
			gen lag = .
			gen number = .
			gen method = "`method'"
			
			local varlist `lagvars' `costlags'
			local listsize: list sizeof varlist
			
			foreach var in `varlist' {
			
				local pos: list posof "`var'" in varlist
				replace var = "`var'" if _n <= (`n'+1)*`pos' & _n > (`n'+1)*(`pos'-1)						
				replace lag = _n - (`n'+1)*(`pos'-1) - 1 if _n <= (`n'+1)*`pos' & _n > (`n'+1)*(`pos'-1)
				replace number = `pos' if _n <= (`n'+1)*`pos' & _n > (`n'+1)*(`pos'-1)
							
				}
			}
		}
			
	end
					
			